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We analyze the quantum phases, correlation functions and edge modes for a class of spin-1/2 and fermionic 
models related to the one-dimensional Ising chain in the presence of a transverse field. These models are the 
Ising chain with anti-ferromagnetic long-range interactions that decay with distance r as 1/r", as well as a 
related class of fermionic Hamiltonians that generalise the Kitaev chain, where both the hopping and pairing 
terms are long-range and their relative strength can be varied. For these models, we provide the phase diagram 
for all exponents a, based on an analysis of the entanglement entropy, the decay of correlation functions, and the 
edge modes in the case of open chains. We demonstrate that violations of the area law can occur for a < 1, while 
connected correlation functions can decay with a hybrid exponential and power-law behaviour, with a power that 
is cr-dependent. Interestingly, for the fermionic models we provide an exact analytical derivation for the decay 
of the correlation functions at every a. Along the critical lines, for all models breaking of conformal symmetry 
is argued at low enough a. For the fermionic models we show that the edge modes, massless for a > 1, can 
acquire a mass for a < 1. The mass of these modes can be tuned by varying the relative strength of the kinetic 
and pairing terms in the Hamiltonian. Interestingly, for the Ising chain a similar edge localization appears for 
the first and second excited states on the paramagnetic side of the phase diagram, where edge modes are not 
expected. We argue that, at least for the fermionic chains, these massive states correspond to the appearance of 
new phases, notably approached via quantum phase transitions without mass gap closure. Finally, we discuss 
the possibility to detect some of these effects in experiments with cold trapped ions. 
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I. INTRODUCTION 

Topological superconductors and insulators have gener¬ 
ated enormous interest in recent years as they correspond 
to examples of novel quantum phases that are not captured 
by the familiar Ginzburg-Landau theory of phase transi¬ 
tions. Breakthrough experiments have already led to the ob¬ 
servation of symmetry protected topological phases both in 
condensed-matter systems [1] and atomic, molecular, and op¬ 
tical physics [2, 3]. While topological phases are finding ap¬ 
plications in fields as diverse as photonics and spintronics, 
the recent probable observation of Majorana modes [4-10] in 
solid-state materials represents the first major step towards the 
realization of topological quantum computing. 

Majorana modes are non dispersive states with zero en¬ 
ergy. In Ref. [11], Kitaev has shown that these modes can 
exist localized at the edges of a one-dimensional supercon¬ 
ductor made of spinless fermions with short-range (SR) p- 
wave pairing. This model is solvable and the underlying lat¬ 
tice Hamiltonian can be mapped exactly onto the well-known 
Ising chain in a transverse field in one dimension. For SR 
interactions, the latter is a text-book example of Hamiltonian 
displaying a quantum phase transition, here from an ordered 
(anti-)ferromagnetic phase to a disordered paramagnetic one. 
Following earlier theoretical works [12, 13], recent experi¬ 
ments with cold trapped ions have generated enormous inter¬ 
est by demonstrating that long-range (LR) Ising-type Hamil¬ 
tonians arise as the effective description for the dynamics of 
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the internal states of cold trapped ions, acting as (pseudo- 
jspins with two or, recently three, internal states. In these 
experiments, effective spin interactions are generated by a 
laser-induced manipulation of the vibrational modes of the ion 
chain [13-17], which are naturally long ranged. The resulting 
Ising-type interactions are antiferromagnetic and decay with 
distance rasa power-law 1 /r"', with an adjustable exponent 
a usually in the range 0 < a < 3.5. 

In experiments with cold ions, the quantum state of individ¬ 
ual particles can be prepared and measured. As a result, both 
the static and dynamical properties of the many-body system 
are accessible. Recent experiments have led to the observa¬ 
tion of instances of interaction-induced frustration [18], non¬ 
local propagation of correlations [19-23] and entanglement in 
a quantum many-body system [16, 24]. Very recently, spec¬ 
troscopy experiments have focused on the precise determina¬ 
tion of the excited states of LR models [25]. 

The experimental works described above are based on the 
understanding of the phase diagram of the Ising-chain in a 
transverse field, which is known exactly for SR interactions 
only. In a seminal work [26], Koffel, Lewenstein and Taglia- 
cozzo have explored the phase diagram of this system with LR 
interactions in the parameter range a > 0.5. The results were 
intriguing: (i) The connected correlation functions decay with 
a power-law tail even within the gapped paramagnetic phase, 
at odds with conventional wisdom inherited from SR mod¬ 
els and consistent with earlier predictions for other quantum 
models with LR interactions [12, 21-23, 27]. Crucially, (ii) 
the entanglement entropy, usually a constant within gapped 
phases, seems to scale logarithmically with the system size 
within the paramagnetic phase for sufficiently small a < 1 as 
well as sub-logarithmically for a > 1. This is remarkable as 
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would signal a violation of the so-called “area law”, dictat¬ 
ing the behaviour of the entanglement entropy in SR quantum 
mechanical systems. These studies also confirm (iii) the per¬ 
sistence of antiferromagnetic and paramagnetic orders with 
decreasing a. 

Research in the area of topological phases with LR in¬ 
teractions is very active, and several possible experimental 
realisations have been recently proposed. In particular, Kitaev 
chains with non-local hopping and pairing may be realized 
in solid state architectures with so-called helical Shiba 
chains, made of magnetic impurities on an s-wave super¬ 
conductor [28, 29]. For atomic and molecular systems, key 
implementations of topological phases have been proposed 
with polar molecules, dipolar ground state atomic quantum 
gases and Rydberg excited atoms [30-47]. In addition, 
the famous Haldane phase may be soon realised in cold 
ion experiments with three internal states per ion [48, 49], 
simulating spin-1 particles. For this latter model, very recent 
theoretical work [50] has demonstrated that major features 
of symmetry-protected topological order can persist for LR 
interactions. It remains an open question to determine the 
validity of these results for generic symmetry-protected 
topological phases with LR interactions. 

In this work, we analyze the quantum phases, correlation 
functions and edge-mode localisation of a class of spin-1/2 
and fermionic models related to the one-dimensional Ising 
chain in the presence of a transverse field. These models 
are the Ising chain with anti-ferromagnetic LR interactions, 
as discussed in Ref. [26], as well as a class of Hamiltonians 
corresponding to a generalization of the Kitaev chain, where 
both the hopping and pairing terms are LR with an algebraic 
decay l/r®, and their relative strength can be varied. 

For these models, we provide the phase diagram for all ex¬ 
ponents a , based on an analysis of: (i) the entanglement en¬ 
tropy; (ii) the decay of correlation functions in all phases; (iii) 
the mass and the localization properties of the edge modes 
when the chains are open. 

In the case of the long-range Ising (LRI) chain we 
utilize numerical calculations based on the density-matrix- 
renormalization group (DMRG) method [51, 52], while 
the long-range Kitaev-type (LRK) models remain exactly 
solvable for all a, allowing for analytical calculation. 

The following results are obtained for all models: 

(i) A violation of the area law for the entanglement entropy 
occurs in gapped regions with a < 1. For a > 1, no violation 
is found. 

(ii) For any finite a, connected correlation functions within 
the gapped phases display a hybrid decay that is exponential 
at short distances and algebraic at long ones. The power of 
the algebraic decay, however, as well as the extension of the 
two decay regimes, depends on a. However, when a < 1, 
the connected correlation functions show a purely algebraic 
decay. 

(iii) For the LRK models, we provide an exact analytical 
expression for the decay of correlation functions within the 
gapped phases that describes the hybrid behaviour with dis¬ 


tance mentioned above and explains its origin. 

(iv) Along the critical lines, we demonstrate that conformal 
symmetry is broken for sufficiently small a, by analyzing the 
finite size scalings of the Von Neumann entropy and of the 
energy density for the ground states, as well as the behaviour 
of the dynamical exponent with a around the minima of the 
spectrum. 

(v) We find the existence of two kinds of edge modes: mass¬ 
less and massive. For the LRK models, massless (Majorana) 
modes, as previously found in [53], appear in the antiferro¬ 
magnetic region of the phase diagram for large a. The antifer¬ 
romagnetic phase for the LRK models is defined in analogy 
with that of the short-range Ising chain. The massive modes, 
instead, are entirely new and are found to appear in a broad 
area of the antiferromagnetic phase for a < 1 and when we 
choose the unbalance e between the strengths of the hopping 
and pairing terms to be different from 1. This results suggests, 
for a < 1, a restoration of the Z 2 symmetry associated with 
the (absence of) ground state degeneracy, as well as a possible 
transition to a novel Z 2 symmetric phase and without mass 
gap closure. Interestingly, if we choose e = 1, the massless 
modes survive for all a > 0 and they are exponentially local¬ 
ized at the edge of the system. 

(vi) For the LRI chain in the antiferromagnetic phase, edge 
modes are massless for all a, and, up to numerical precision 
are exponentially localized at the edge of the chain, in con¬ 
trast, e.g. to [50]. However, surprisingly we find in the param¬ 
agnetic phase a localization of excited, gapped, energy eigen¬ 
states for a < 1, which for a > 1 are instead delocalized in 
the bulk. 

(vii) We finally discuss the persistence of some of the LR 
effects discussed above (e.g., hybrid decay of correlation 
functions and edge mode localization) in small chains of up 
to 30 sites, as relevant to current experiments. 

The paper is organized as follows. In Sec. II we intro¬ 
duce the model Hamiltonians that we consider in this work 
[Sec. II A], the observables that are used to characterize the 
various phases [Sec. II B], and present the corresponding 
phase diagrams [Sec. II C]. In particular, in Sec. IID we dis¬ 
cuss the critical lines of the Ising and Kitaev models, and ar¬ 
gue that conformal symmetry is broken for sufficiently small 
a. In Sec. Ill we provide an analytic calculation of the corre¬ 
lation functions for the fermionic Hamiltonians that explains 
the hybrid exponential and algebraic decay observed in these 
LR models [Secs. IIIA and IIIB], In Sec. HID we provide 
a numerical comparison with results for the LRI chain, dis¬ 
playing similar behaviour. In Sec. IV we analyze the edge 
modes in the LRI and LRK chains. In particular, in Sec. IV A 
we analyze the properties of gapless Majorana modes that are 
found in the anti-ferromagnetic phases of the LRK models for 
a > 1. In Sec. IV B, instead, we demonstrate that the edge 
modes can become massive for a < 1, signalling a transition 
to a new phase. We discuss similar results obtained for some 
excited states that get localized on the edges in the paramag¬ 
netic phase of the LRI model for a < 1. In Sec. V we discuss 
the observability of some of the results above in small chains 
of a few tens of particles, as relevant for cold ions experi- 
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ments. Finally, Section VI discusses conclusions and outlook. 


II. MODEL HAMILTONIANS AND QUANTUM PHASES 

In this section we introduce the model Hamiltonians that 
we consider in this work and present the corresponding phase 
diagrams that we compute based on results from the entan¬ 
glement entropy, decay of correlation functions, spectrum of 
excitations, and edge mode localization, as discussed in detail 
in the following sections. 


A. Model Hamiltonians 

1. LR Ising model. 

In this work, we are interested in Ising-type Hamiltonians 
with LR interacting terms. The long-range Ising Hamilto¬ 
nian [26] reads 


crlcrl J' ' 

// LR i = sin0 ^ ' ! a + cos8^a~, (1) 

i=l,j>i ' l i=l 

where cr v . (v = x, y, z) are Pauli matrices for a spin-1/2 at site j 
on a chain of length L. The first term on the right hand side of 
Eq. (1) describes spin-spin interactions that we choose antifer¬ 
romagnetic (AM) with sinf? > 0 (or equivalently 0 < 8 < n). 
The second term describes the coupling of individual spins to 
an external field pointing in the z-direction. Thus, while the 
first term favors an antiferromagnetic phase with spins point¬ 
ing along the x direction, the second terms favors a paramag¬ 
netic (PM) phase where all spins align along z. In the case of 
SR interactions (i.e., for a — > oo) the Hamiltonian Eq. (1) is 
exactly solvable and a quantum phase transition between these 
two phases is known to occur at 8 C = n/4. Reference [26] has 
shown numerically that a quantum phase transition separating 
the AM and the PM phases survives for all finite a > 0.5. 
Below, we are interested in exploring the phase diagram of 
Eq. (1) for all a and 8. 


2. Long-range Kitaev chains. 


Related to the LRI chain, in the following we introduce and 
analyze a class of fermionic Hamiltonians of the form 


Hlrk = sin 8^ 


aj aj + (1 + e) ciiCij + h.c. 


+ 2 cos 




( 2 ) 

where a t describes the creation operator for a fermionic par¬ 
ticle at site j, and tij = ajaj. The Hamiltonians Eq. (2) repre¬ 
sent generalizations of the Kitaev chain for spinless fermions 
with superconducting p-wave pairing, where both the hopping 
and pairing terms decay with distance algebraically with ex¬ 
ponent a. Here, all energies are expressed in dimensionless 


units and the parameter e governs the unbalance between the 
hopping and pairing terms. 

In the SR limit a —» oo, Eq. (2) maps into the Ising chain 
Eq. (2) via the Iordan-Wigner transformation [54] 


o-j = a . exp 


V (<j ) 


cr) = 2n, - 1, 

J J 


(3) 

(4) 


with crt = (<xj+icrp/2. However, at finite a this identification 
does not hold anymore due to the contributions of the string 
operators exp(i 7tYit<jne) in Eq- (4). In particular, unlike the 
LRI, Eq. (2) remain exactly solvable, allowing for analytic 
solutions at any finite a. 


B. Observables 

1. von Neumann Entropy. 

Entanglement measures are routinely used to characterize 
the critical properties of strongly correlated quantum many- 
body systems [55]. A key example is the von Neumann en¬ 
tropy S( that we employ in this work. For a system of L sites 
that is partitioned into two disjoint subsystems A and B con¬ 
taining t and L — € sites, respectively, Sf is defined as 


S ( = -Trpf log 2 pt. 


(5) 


where pi is the reduced density matrix of the subsystem A. 

Two general behaviors of Sf are known for the ground 
states of one-dimensional SR interacting systems. Within 
gapped phases, Sf saturates to a constant value independently 
of ( and thus obeys to the so-called area law [56], On the 
contrary, Sf diverges with f for critical gapless phases and, 
for conformally invariant systems, satisfies the universal law 
[57]: 


S t = a + log, d({, L), (6) 

6 

with d({,L) = (2L/n) sin(n{/L) for the case of open bound¬ 
aries. Here, a is a non-universal constant and c is the cen¬ 
tral charge of the theory. The latter characterizes the univer¬ 
sality class of the gapless phase [58, 59]. Sf, and thus the 
central charge, can in principle be directly computed by nu¬ 
merical techniques, such as Density Matrix Renormalization 
Group [51, 52] as well as by means of analytical methods for 
quadratic Hamiltonians [60-63], 

In the case of LR models it has been shown [26, 53, 56] 
that the divergence of Sf in Eq. (6) can also occur for gapped 
phases, corresponding to a so-called violation of the area law. 
Since this violation is found to be logarithmic, an effective 
central charge c e ff may be defined also within the gapped 
phases and used to characterize the main features of the phase 
diagram for LR interactions [26], 
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2. Correlation functions. 

The various quantum phases can be characterized by the 
decay of two-point correlation functions with distance. For 
the LRI model, we are interested in the connected correlations 

C% = (o y i o!)-(o?)((r]), (7) 

with v = x, z. 

For the LRK models we are interested both in the density- 
density correlation function 

gi(i, j ) = <«;» j ) - i>h) (n j) ( 8 ) 

and the function 


p- wave SR pairing interactions [described e.g. by the limit 
a —> oo of Eq. (2)]. The existence of these localized states is 
related with the spontaneous breaking of the discrete Z 2 sym¬ 
metry associated with the parity of the fermion number: when 
this Z 2 symmetry is broken, two degenerate ground states 
with different parity appear. Here, they will be labeled by |0 + ) 
and |0 } in the even and odd parity sectors, respectively. 

As Hamiltonians Eq. (1) and (2) are equivalent in the limit 
a —> oo, the same breaking of Z 2 symmetry (now related to 
spin flips along x direction) described above occurs also for 
the LRI chain, resulting in the presence of two degenerate 
edge states in this model. From the discussion above it turns 
our that the analysis of the edge modes can be utilized to 
characterize the quantum phases of the system. 


£(/, j) = ((a] + a t ) exp 


7T1 


j ~ l 

i 2 a\a t 


I (=i +1 


(a) + aj)) 


(9) 


that correspond to the functions C zz and ( cr x cr x ) in the LRI 
model, respectively, via the Jordan-Wigner transformation 
given above [54], Since the LRK models are quadratic, the 
functions Eq. (8) and (9) can be directly obtained from the 
one-point correlations ( a\aj } and (a\a\) via Wick’s theorem. 
In particular, one finds 


£ 2 ( 1 , 7 ) = \(a]a])\ 2 - \(a]aj )| 2 
and 


S(/, j) = det 


Gi,i +1 
G 7 -U+I 


Gij 


G i-hj 


with G m ji = 6„ 


+ 2 (a m a n ) + 2 ( a m a n ). 


( 10 ) 


( 11 ) 


For SR interactions, the connected correlations above are 
known to decay exponentially (algebraically) with distance 
within the gapped (gapless) phases. Surprisingly, for several 
models with LR interactions it was reported that algebraic 
decay of correlations can coexist with an initial exponential 
decay within gapped phases [12, 22, 53], An analytic 
understanding of this effect has so far proven elusive. 

In the following we use the decay of correlation functions 
to characterize the various phases. In particular, for the LRI 
chain we provide extensive numerical results using DMRG 
techniques, while for the LRK chains we exploit the integra- 
bility of the models to derive an analytic expression of the 
behaviour of correlation functions in all parameter regimes. 


3. Edge states and edge gaps. 

Localized edge states within topological phases have at¬ 
tracted much interest over the last decade, largely because 
of possible applications in schemes for topological quantum 
computing [64], In particular. Ref. [11] has shown that lo¬ 
calized states arise at the edge of a ID superconductor with 


In this work, we identify the localized states by directly 
computing their wavefunctions and masses, which can be 
achieved either numerically for Eq. (1) using DMRG simu¬ 
lations or exactly for Eq. (2). In order to accomplish this 
task for the LRI model, it is useful to exploit the Jordan- 
Wigner transformation Eq. (4) to define new fermionic op¬ 
erators cj = cr + j exp(i7iY J c<i crJ i ar ~{) [similar to Eq. (4)] from 
spin operators cr*. We then compute the wave-functions ip^ ' 2) 
of the massless edge modes as [65] 

<pf = (0 _ |cj|0 + ) + <0 + |c}|(T> (12) 

^ 2) =i((0-|c]|0 + )-<0 + |c]|0->). (13) 

Here, the states |0 + ) and |0~) are the ground states in the even 
and odd parity sector also for the LRI model, respectively. 
The mass of the two modes tp. , also known as edge gap (in 
order to distinguish it from the usual mass gap, i.e. the energy 
difference between the ground state and the first excited 
bulk state, see Refs. [50, 66]), is defined as the difference 
Flo-) - £|o+>, with £|o±> being the energy of the state lO*). 

In the following we will also be interested in characterizing 
the localisation of massive edge states that are found in the 
phase of LRI where the Z 2 symmetry is preserved, and thus 
where a unique ground state |0 + ) occurs in the even-parity sec¬ 
tor. This localization arises for the first two excited states |1 _ ) 
and 2 ) for a < 1, which are instead delocalised in the bulk 
for a » 1 [see Sec. IV B], Their wavefunctions read 

wf - (r|c]|0 + > wf - (2~|cj|0 + ). (14) 

Here the mass of the mode wj 1 is defined as the difference 
A E(s) = Lij-) - £|o+) with s = 1,2. The difference between 
Eqs. (12), (13) and (14) is that in the first two expressions the 
second term on the right hand side of the equations is nonzero 
because of the zero-energy condition [65]. 

For the LRK models, the wavefunctions for both the mass¬ 
less and massive modes can be extracted following Ref. [54]. 
The latter describes a technique for the exact diagonalization 
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of a generic fermionic quadratic Hamiltonian of the form 




U= 1 


Aijalaj + 2 


^ («; Bija] + H.c.) 


(15) 


with A = [Ay] and B = [By] real matrices. Equation (15) can 
be cast in diagonal form as 


L-l 

H = ^ j A n T]lTj n . ( 16 ) 

n =0 

by a singular value decomposition of the matrix A + B: 

A„ = ^ j(Aji + B'ji)<p in . (17) 

<j 

Here, A„ are single-particle energies (ordered as Ao < Ai < 
• • • < A | J and q,, are fermionic operators defined by the 
following Bogoliubov transformation 

q n = ^ g n ja j + K ja). (18) 

i 

The matrix elements <p„j = g„j + h„j and \jt n j = g n j - h n j can be 
directly identified as the wavefunctions of the two Majorana 
modes ( 77 ,, + qj,)/2 and i(q n - qj,)/2 with energy A,„ while g nj 
and h n j are the wavefunctions of q„ and q n . 


C. Phase diagrams of LRI and LRK models 

In this section, we present the phase diagrams for the LR 
models Eqs. (1) and (2), obtained from an analysis of the ob¬ 
servables described above. The results for the LRI model were 
obtained numerically via DMRG techniques for chains of a 
length L up to L — 200. For all calculations we utilized up 
to 128 local basis states and 10 finite-size sweeps [51, 52], 
The discarded error on the sum of the eigenvalues of the re¬ 
duced density matrix was always less than 10 8 . For the LRK 
models all results were obtained (semi-)analytically. 


1. LR Ising model. 

Our results for the phase diagram of the LRI model are 
summarized in Fig. 1(a), where we plot the effective central 
charge c e ff defined in Sec. IIB 1 as a function of the angle 
8 and the power a of the antiferromagnetic term in Eq. (1). 
Hamiltonian Eq. (1) is invariant under the transformation 
cos 8 —> -cos 8 and thus the phase diagram is symmetric 
around 8 = n/2. 


find that the behaviour of the full correlation functions (cr x a x ) 
and (trier:) is consistent with the persistence of paramagnetic 
and antiferromagnetic orders for all a. However, different 
from the SR model, we find that the connected correlation 
functions decay with distance with a hybrid behaviour that 
is exponential at short distances and algebraic at long ones. 
An example is shown in Fig. 2(a) for C XX R in the PM1 phase. 
Surprisingly, we find numerically that the exponent y x of 
the long-distance decay for C XX R displays three difference 
behaviours: (i) for a > 2 it fulfills y x = a, consistent with 
the results of Refs. [12, 26]. However, (ii) for 1 < a < 2 
we obtain a hybrid exponential and algebraic decay with a 
different y that depends linearly on a with a slope consistent 
with ~ 0.55(5) and (iii) for a < 1 we observe numerically 
a curve compatible with a pure algebraic decay, with an 
or-dependence of y x that is linear with slope ~ 0.25(2). The 
fitted exponent y x is shown in Fig. 2(c). 

For a < 1 in the paramagnetic regions of the phase diagram 
denoted as PM2 we find that the effective central charge grows 
continuously with decreasing a from zero to a finite value that 
appears to be 0 -dependent and has a maximum of order 1 for 
a = 0 and 8 ~ n/2. An example for 8 = 0.2n is plotted in 
Fig. 1(d) (blue squares). As mentioned above, in this PM2 
region, the correlation function C XX R is found numerically to 
decay as an almost pure power-law. 

The energy spectrum changes in this region PM2 compared 
to the case PM1: the energy gaps AF(1) and AE(2) between 
the ground state and the first excited states in the odd parity 
sector increase with decreasing a , as shown in Fig. 8 (b), and 
the wave-functions of the two lowest-energy excited states 
|1 _ ) and |2~), defined in Sec. IIB 3, become localized at the 
edges of the chain [Fig. 7(b)]. 

In the antiferromagnetic phase, denoted as AM, the effec¬ 
tive central charge is instead zero for all a [26], For a > 1 the 
connected correlation functions Cf R display a clear algebraic 
decay at long-distances [see the example in Fig. 6 (a) below], 
while our numerical results do not allow for establishing 
whether an initial exponential decay is also present, as 
expected. The ground-state is found to be doubly degenerate 
for all a. This degeneracy is due to the spontaneous breaking 
of the Z 2 spin-flip symmetry [58, 67], and is related to the 
two modes that are localized at the edges of the chain as in the 
short-range Ising model [ 68 ]. While a LR power-law tail may 
be present [50], the localization of these modes is here found 
to be consistent with exponential up to numerical precision 
[see Fig. 7(b) below]. We come back to this point below. 


2. LR Kitaev models. 


We find that for a > 1, t: c ii is zero everywhere except along 
two critical lines. By comparing with results for the energy 
gap (not shown), we find that the critical lines separate two 
gapped regions [denoted as PM1 and AM in Fig. 1(a)] that for 
a —> 00 correspond to the known paramagnetic and antiferro¬ 
magnetic phases of the SR model. Similar to Ref. [26], we 


The phase diagram of the LRK model Eq. (2) for e = 0 
is reported in Fig. 1(b), where we plot the effective central 
charge c e ff defined in Sec. IIB 1 as a function of the angle 
8 and the power a of the decay of the pairing term. In this 
case the invariance of Eq. (2) under cos 8 —> - cos 8 is lost for 
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FIG. 1. (a) Phase diagram of the LRI model Eq. (1) from the effective central charge c e fr extracted from the entanglement entropy Si [see 
Eq. (6)] for a system of L = 100 sites, as a function of a and 9. The phase diagram is symmetric with respect to 8 = nj2. PM1 and AM denote 
gapped paramagnetic and antiferromagnetic phases, respectively, with c e fr = 0. In the paramagnetic gapped region PM2, c e ff > 0 [see panel 
(c) blue squares] and massive edge modes [see Sec. IV B and Fig. 7(d)] are found. The continuous black lines correspond to critical transition 
lines with central charge c = 1/2. The dashed red lines signal an increase of the central charge from 1/2 up to c ^ 1 with decreasing a [see 
Fig. 3(b) blue triangles], (b) Phase diagram of the LRK model Eq. (2) from c e ff extracted from St in the thermodynamic limit as a function of 
a and 8, and for e = 0. AMI and PM are gapped phases with c cff = 0 for a > 1 and c e ff t 0 for a < 1 [see panel (c) green circles], AMI shows 
massless edge modes with a hybrid exponential and power-law decay with distance. AM2 shows edge states with a finite mass Ao, which 
increases with decreasing a [see panel (c) and Fig. 8(a)], The continuous black lines correspond to critical transition lines with central charge 
c = 1/2. The dashed red lines signal an increase of the central charge from 1/2 up to c - 1 with decreasing a. (c) Effective central charge c e ff 
for the LR Ising (blue squares) and the LRK models (green circles) for 9 = 0.2n as function of a. (d) Scaling with the system size L of the 
mass A 0 of the edge modes within the AM2 phase, for 8 = 0.15n and for several a < 1. (e) Edge gap A 0 of the LRK model for a system of 
L = 400 sites and e = 10. In the region denoted by “edge gap" the Zi symmetry of the model is restored and the AM2 phase appears. 


any finite a and the phase diagram is not symmetric around 
8 = 7t/2. 

Figure 1(b) shows that for a > 1 and 0 < 8 < n, two phases 
exist that are denoted as PM and AMI separated by two 
critical lines. In the limit a —* oo these phases correspond to 
the paramagnetic and antiferromagnetic phases of the LRK 
model and are gapped. Consistently, Fig. 1(b) shows that the 
effective central charge c e ff is zero within the phases for all 
a > 1 [see Fig. 1(d) for an example], while c e ff + 0 along 
the critical lines, as expected from general results for SR 
systems [57]. 

The PM and AMI phases are distinguished by different 
asymptotic values of the correlation functions 'L(i, j) defined 
in Sec. IIB2. In the region denoted as AMI, Z(i, j) has a 
finite value for | i - j\ —> oo, while £(/, j) decays for | i — j\ —> oo 
within the PM phase. Similar to the situation in the LRI 
model (see above), the decay to zero of X(i, j) in the PM 
phase shows a hybrid exponential and power-law behaviour 
with distance. This is shown for a particular value of 8 in 
Fig. 2(b), where we find numerically that the exponent for 
the power-law tail of X(i, j) equals y z = a when a > 1. For 
a < 1, however, the exponential part becomes numerically 
unobservable, and £(/,_/') decays essentially algebraically 
within the PM phase with an exponent that grows to 2 for 
a —> 0 [see Fig. 2(b)]. 

Remarkably, we show below in Sec. Ill that the hybrid ex¬ 
ponential and algebraic behaviour described above can be ob¬ 
tained analytically in all phases for several correlation func¬ 


tions, such as the one-body and the density-density correla¬ 
tion functions. In particular, the leading contribution to the 
one-body correlation function (aj aj ) reads 


(a R ao) = m a ,e ■ -p-1- Sff.e ■' 


1 


R a+ 1 

1 

R 2a -' 

1 


ri¬ 


ch > 2 

l<a<2 (19) 
0 < a < 1 


where the pre-factors /?(„// and 'B,,# are derived below in 
Sec. III. The algebraic part of the decay of g 2 (R) is instead 
found to be gn(R) ~ R la and giiR) ~ R 2 for a > 1 and 
a < 1, respectively. 


In Sec. IV A, we show that a similar hybrid exponential 
and power-law decay is found for the localization of the edge 
modes within the antiferromagnetic phase AMI: For a > 1, 
the edge modes are massless, as expected from the SR Kitaev 
model. However, for a < 1 the edge modes acquire a finite 
mass, i.e., become gapped. This is shown in Fig. 1(c), where 
we plot the edge gap Aq defined in Sec IIB 3 as a function of 
a for a few values of the parameter e of Hamiltonian Eq. (2): 
while for a > 1 the gap scales to zero with the system size 
L as 1 /L“, for a < 1 it remains finite and for a ~ 0 can be 
of order unity. The presence of the gap removes the degener¬ 
acy of the ground-state, signaling a new phase for this class of 
topological LR models. This latter phase is denoted as AM2 
in Fig. 1. 
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with 

A a (k) = 2 sin 8 ^(g a (k) + cot 9) 2 + (\ + e) 2 f a (k) 2 . (22) 

Here, the new fermionic operators (rf k , ij k ) are given in terms 
of the operators (a k , a _*) by 

4 ) = ( sin & icosyJtW a\ \ 
n-k) \!COS p k sin p k j\a- k J’ 

where tan (2/^) = (1 + e)f(k)/[g(k) + cot 6], The ground 
state of Eq. (21) is the vacuum of the rj k fermions and has an 
energy density e 0 (a, L) = - % k A a (k)/(2L). 

The critical lines can be computed from the dispersion re¬ 
lation Eq. (22) as follows: (i) For a finite system with L sites, 
Eq. (22) is zero on the line cot 6 + g a (0) = 0 and on the line 
cot 0 +# 0 ,( 71 ) = 0; (ii) For a system in the thermodynamic limit 
L —» 00 instead f a (k) = Tm Li„(e 1<; ) and g„(k) = Re Li„(e l#; ), 
Lio(x) being the polylogarithm of order a [69], and the critical 
line cot 6 + g a ( 0 ) = 0 ends at the point 6 — n and a = 1 . 


FIG. 2. (a) C“ correlation Eq. (7) for Hi sing (0 = 0.2tt and L = 
60), showing the hybrid exponential and power-law behavior for a > 
1 and the purely power-law for a < 1. The oscillations and the 
bending of CA, for R > 50 are due to finite size effects of the DMRG 
simulations, (b) 1.(1, R) correlation Eq. (9) for // LRK at 9 = 0.2 n 
and L = 100 (PM phase), showing the same hybrid behavior and 
the same decaying exponent as C* X R for a > 1. For a < 1 instead 
the power-law tails have exponents differing from C\ X R . (c) Decaying 
exponent y x of the algebraic tail of the C\ X R correlation fitted as 1 /R r -'. 
Three different behaviours for y x , corresponding to the three dashed 
black lines are numerically observed (see Sec. IIC 1 ). (d) Decaying 
exponent of the algebraic tail of the E(l.R) correlation fitted as 
1 IF™. 


D. Critical lines 

1. LR Kitaev models. 


For the LRK models, the critical lines can be computed ex¬ 
actly as Eq. (2) are integrable. A Fourier transform of the 
fermionic operators a t takes Eq. (2) to the form 


H L rk = sin0 V (, 


k 


&k Cl-k 


)({ 


ga(k) + cot 6 
i(l + e)f a (k) 


-i(l + e)f a (k) \ / a k 


-(g a (k) + cot 0)] \a 


-k. 


( 2 °) 

where a k reads a k = 2 / c' k 'cij, k = Imi/L is the lattice 

momentum with n = 0,..., L — 1, and the functions f„(k) and 
ga(k) read f a (k) = £ f sin (k£)/£ a and g a (k) = Xe cos(k£)/{ a , 
respectively. 


A Bogoliubov transformation brings Eq. (20) in diagonal 
form as 


//|,rk - a (k) (4m - 1 /2), (21) 

k 


For the critical line with 9 < zr/2 we compute the value 
of the central charge c by two methods: (i) by fitting the von 
Neumann entropy Eq. ( 6 ) and (ii) by studying the finite size 
corrections to the ground state energy density eo(a,L) (see 
Ref. [53] for a similar model). 

The results for c obtained from the scaling of the von Neu¬ 
mann entropy Eq. ( 6 ) are reported in Fig. 3(b). For a > 1, we 
find c = 1/2 as expected from the SR model. For a < 1, how¬ 
ever, c e ff increases up to values of order one [see red dashed 
line in Fig. 1(a)]. In Ref. [53] it was demonstrated for a related 
model with LR pairing only that this behaviour corresponds to 
an exotic change for the decay of density-density correlation 
functions: For a < 1 their oscillations mimic those of a Lut- 
tinger liquid. Here, we find a similar behaviour (not shown). 
The increasing of c e ff, below a - 1, is also found in the very 
recent work [70] where the scaling of the von Neumann en¬ 
tropy in the thermodynamic limit is analytically analized. 

This anomalous behaviour of c e fr points towards a breaking 
of conformal symmetry along the critical line, which we 
analyze further below. 

The breaking of conformal symmetry can be inferred also 
by analyzing the scaling of the energy density eo (a, L) with 
the system size L [53]: For a conformally invariant theory the 
following relation must hold [59] 


7r vf c 

e 0 (a, L) = e 00 (a) - 2 , (24) 

where e C o(o') is the energy density in the thermodynamic 
limit, Vf = |2 sin 0Li,,_ 1 (-1 )| is the Fermi velocity and c is the 
central charge of the conformal theory. We analyzed numeri¬ 
cally eo (a, L), finding that relation (24) works properly only 
for sufficiently large a > 2. This results in a value c - 1 /2 
for a > 2, as expected from results for the SR Kitaev chain. 
Conversely, for a < 2, eo (a, L ) does not satisfy the scaling 
law (24), which implies a breaking of conformal symmetry. 
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This behaviour also implies that the quantum phase transition 
between the PM and the AMI phase for 8 < n /2 and a <2 
is in a different universality class from that of the SR Kitaev 
(Ising) model. We notice that, even if conformal symmetry 
is broken, the von Neumann entropy predicts a value for 
c e ff which tends to 1 as a goes to zero, compatible with the 
observed decay of density-density correlation functions. 

We further confirm the breaking of conformal symmetry for 
the fermionic models by looking at the behaviour of their low- 
energy spectra. The dispersion relation of a conformal field 
theory is linear in the momentum k , implying a dynamical ex¬ 
ponent z = 1 [59], Consistently, by expanding the dispersion 
relation A a (k) for the long range Kitaev models on the critical 
line with 8 > n/2, for a > 2 we find A a (k ) ~ k, for k —» 0. 
However, for 1 < a < 2 we obtain the scaling A a (k ) ~ k' r 1 . 
This latter scaling implies a dynamical exponent z = a — 1 
that varies continuously with a and is different from that of 
a conformal held theory. This would imply that the quan¬ 
tum phase transition between the PM and the AMI phase for 
8 > 7t/2 and a < 2 is in a different universality class from that 
of the short-range Kitaev model. The appearance of a new 
universality class due to long-range interactions is also found 
in Refs. [71, 72]. Incidentally, we notice that linearity of the 
spectrum around the minimum is only a necessary condition 
for the persistence of conformal invariance: indeed along the 
critical line at 8 < zr/2 conformal symmetry breaking arises 
even if the low-energy spectrum is linear for every a. 


2. LR Ising model. 

We locate the critical line (for 8 < 7r/2, the other being 
symmetric) of the LRI model numerically by using two 
complementary ways that agree up to finite-size effects. 
Firstly, we determine the points in the phase diagram of 
Fig. 1 (a) where the energy gap between the ground state and 
the first excited state reaches its minimum. Then, we compute 
the effective central charge c e ff for different system sizes L 
from Eq. (6) and determine the precise values of a and 8 for 
which its value does not depend on L [73, 74], Examples of 
this latter technique, which is found to be particular precise, 
are presented in Fig. 3 for different a. We notice however that 
this method does not allow us to extract a precise value for 
c e ff when a < 0.2, since within our numerical results lines 
with different L do not cross at a single point in this region. 

On the critical line, we find that for a > 1 (black solid 
lines in Fig. 1) c e ff is equal to 1/2 as expected for the central 
charge of the critical SR Ising model. However, for a < 1 (red 
dashed lines in Fig. 1), c increases continuously up to a value 
of order 1 as shown in Fig. 3(b). We argue that on this line 
the conformal invariance of the model is broken as the found 
values of c do not coincide with the discrete set allowed for 
the known conformal field theories [58, 67]. 

Based on the mismatch between the predictions for c from 
the von Neumann entropy and the ground state energy den¬ 
sity found in the previous subsection for the LRK model, we 



FIG. 3. (a) Effective central charge extracted from the scaling with 
L of the von Neumann entropy (6) as function of 6 for two values of 
a, plotted for different system sizes L = 60, 80, 100. We located the 
critical point (indicated by red triangles) where the effective central 
charge does not depend on L. (Main panel) scaling at a = 2.0: the 
critical point is 6 C ~ 0.3027T and the central charge c is compatible 
with 1/2. (Inset) scaling at a = 0.7: the critical point is 8 C » 0.3767T 
and the central charge c is compatible with c ~ 0.66. (b) Central 
charge on the critical lines of the LRI model (blue triangles) obtained 
with the same method as in panel (a) and of the LRK model (green 
diamonds) obtained by finite size scaling from the expression for en¬ 
tanglement entropy in Eq. (6). 



FIG. 4. (a) Integration contour adopted to evaluate the integral (25). 
The red dashed line is the branch cut of the square root in the denom¬ 
inator of the integrand in (25). (b) Decay rate 1/^ of the exponential 
part of the correlation of Eq. (29) as function of a for two values of 
6. The red lines represent the value of l/£ for a —> oo. 


cannot exclude here a conformal symmetry breaking also in a 
certain range for a above a = 1. However our results do not 
allow us to provide a final answer, since our DMRG calcula¬ 
tions for the LRI model cannot reach sufficiently large sizes to 
perform a satisfying finite-size scaling for the energy density. 


III. CORRELATION FUNCTIONS FOR THE LRK 


In this Section, we present an analytical calculation of the 
one-body correlation functions for the LRK models. The 
latter display a hybrid exponential and algebraic decay with 
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distance that is explained by exploiting the integrability of 
the models. Higher-order correlation functions, such as the 
density-density correlations, are readily obtained from these 
correlations via Wick’s theorem [see Sec. IIC 2 and below]. 

The one-body correlation functions (a R ao) and (ajjal) read 



respectively. In the following, we focus first on the one-body 
correlation function (a R a ( )} and come back to the anomalous 
correlation (a R a (j ) in Sec. IIIC below. 

In order to evaluate Eq. (25), we use the Cauchy theorem 
applied to the contour in the complex plane drawn in Fig. 4 

< “ ;oo>= -h (i + L + L + IJ e “ (27) 

with Q a (z) = [ga(z) + cot#]/[2 7(cot # + g a (z)) 2 + f 2 (z)] and 
z = k + iy, and where we have chosen e = 0. In Eq. (27) 
we have neglected the contributions from C± and C'_ as they 
vanish when M —> oo. As we explain below, the integrations 
over the lines £ ± and Co, 2 n (and thus momenta k - n and k - 0 
and 2n, respectively) are responsible for the exponential and 
algebraic behaviour observed in these models, respectively. 


A. Exponential decay 


The sum of the integrals on the lines (where z - n + iy) 
and £+ (where z - n + + iy) of Fig. 4 gives 

gi nR g -£R 

1? P (R) = - e-> R g a (n + \(y + {))dy. (28) 

n Jo 

Equation (28) displays an exponential behaviour with a decay 
constant 1 j£. The appearance of this quantity is due to the 
square root in the denominator of @ a (z), yielding a branch cut 
from z,\ = 7r+i £ to oo. The leading term in Eq. (28) is obtained 
by integrating @ a (7t+i(y+£;)) in the limit y —> 0 [75], and reads 

e -f« 

/“ p (R) = (-l) s ^, e —, (29) 

with 

(2 cot # + Li,/-e f ) + Li,/-e“ f )) 

J\ a o =-77T- . (JU) 

4 Vtf|Lia-i(-ef)| v [cot6> + Li^-e-f)] 1 ^ 


The decay constant is related to the zeroes of the denom¬ 
inator of @ a (z) and is obtained by solving the equation 


[cot# + Li^-e ^)][cot# + Lio(-e^)] = 0. 


Two cases must be distinguished: if cot# < |Li c /— 1)[, the 
equation [cot# + Li^-e - ^)] = 0 admits a solution, since 
the function -Li n (-e - ' f ) is always decreasing for £ > 0. If 
instead cot# > |Lio(—1)|, then [cot# + Li^-e^)] = 0 admits a 
solution, for the same reason as above. In the following we 
focus, without loss of generality, on the first case, where £ is 
solution of cot# = -Li ( /-e - ^). 

Notably for a — 0 solutions exist only for # > 7r/4 [in which 
case, f = +log(tan# - 1) for cot# § 1/2]. For # < tt/4, 
instead, Eq. (31) does not admit any solutions, which implies 
the absence of exponential decay. This is in contrast with, 
e.g., the expected behaviour of correlation functions within 
gapped phases for SR models. 

Figure 4(b) shows the decay constant 1 /f as a function of 
a for two different values of #. In particular, for # = 0.27T 
and (r —» oo, l/£ tends to the SR value (£ —» |log|tan#||), 
as expected. However, for a —> 0 we find that 1 j£ tends to 
zero, essentially linearly with a. As explained below, this can 
result in the non observability of the exponential dependence 
of correlation functions for a < 1. Notice however that even 
if £; is finite when a < 1 for # = 0.47T, the exponential part of 
the correlation functions is still unobservable. 


B. Algebraic decay 


The sum of the integrals on the lines Co (where z = # + ry) 
and Con (where z = 2n - r\ + iy) of Fig. 4 gives 

1 f°° 

ir(R) = - e^ s Im^(iy)dy, (32) 

n Jo 

after sending q 0. The leading contribution for R —> oo to 
the integral in Eq. (32) can be computed again by integrating 
the imaginary part Im Q„(\y) for y —> 0. By exploiting the 
following series expansion of the polylogarithm [69] 


Lp/e^') = T(1 - a)(+yf~' 



£(a ~ j ) 
./! 


( ±y) ] , 


(33) 


one obtains 


with 


7 pow (R) = B aJ!) • - 


R a+1 

a > 2 

1 

1 < a < 2 

Rla-\ 

1 

R 2 ~ a 

a < 1 


(34) 


a£(a - 1) 

- w a > 2 

(£(a) + cot #) 2 

(cos(Tra) - 1)T(1 - a)F(2a - 1) , 

--- 1 < a < 2 

4T(a)(cot(#) + ((a)) 2 

cos 3 (f)T(2 - or)(cot(#) + ((a)) 


(31) 


7rT(l - a) 


a < 1 
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FIG. 5. (a) Comparison between numerical solutions of Eq. (25) 

(symbols) and the analytic expressions Eqs. (29) and (34) (solid 
lines) for 6 = 0.239^ and different a. The inset refers to a = 0.45, 
where the exponential part is negligible, (b) Cf +1 R ^ +R correlation for 
Hlri at 9 = 0.207n, L = 100, R 0 = L/4 and different a. The power- 
law part of Cj? +1R +R shows the same decay exponent y z = 2a for 
a > 1, as g2(0, R) for Wlrk (Eq. (39), solid black lines). 


correlation function gi(i, j ) reads gi(i, j ) = [(a. a j.)| 2 - |(a ( ! aj)\ 2 
and, from Eqs. (34) and (36), the leading part of its LR power 
law tail is found to be 


£2(0, R) 


-z- a > 1 

R 2a 


~ R 2 0<a<l. 


(39) 


We notice that a hybrid exponential and algebraic behavior 
similar to that described above has been already observed 
numerically in certain spin and fermionic models, e.g., in 
Refs. [12, 26, 50, 53], This behaviour is characteristic of 
the non-local interactions and from our analysis appears to 
be largely unrelated to the presence or absence of a gap in 
the spectrum. In fact, we have shown here for Hamiltoni¬ 
ans Eq. (2) that this hybrid decay does not require exotic 
properties of the spectrum (see Refs. [21] and [22]), rather 
is due to the different contributions of momenta k = n (as 
for the SR limit) and k = 0,2zr respectively. In particular, 
the latter momenta are responsible for the LR algebraic decay. 


While, e.g., the phase diagram of Fig. 1(b) demonstrates the 
persistence of individual paramagnetic and antiferromagnetic 
phases with varying a, the analytic expressions Eqs. (34) 
and (35) for the one-body correlation function clearly show 
that different regions are in fact present within each phase. 


C. Hybrid decay and other correlations 


The two contributions from Eqs. (29) and (34) sum up to 
give the hybrid exponential-algebraic behaviour of Eq. (19) 
valid for all a. Figure 5(a) shows that the analytical results 
are in perfect agreement with a numerical solution of Eq. (25). 


The anomalous correlation (u R a (} ) can be computed along 
the same lines as before and is given by 


( a R a q) - 3\ q 


,(-l) 




Vr 


+ £' 


>( a ) 


— a > 1 
R a 


- 0 < a < 1 

R 


(36) 


with 


3K 


(«) _ 


Lii-e-f) - Lii-ef) 


a ’ e 4 |Li a _ 1 (-ef)| 1/2 (cot6 | + Lii-e~f)) 1/2 


(37) 


and 


S (a) = 


1 


4(cot0 + {(a)) 


a > 1 


lna 

-cos — 0 < a < 1. 

2n 2 


(38) 


From (a^ao) and W R a o) we can compute other local corre¬ 
lations by Wick’s theorem. For instance, the density-density 


Finally, as mentioned before, we notice that (i) the contribu¬ 
tion to the imaginary part in Eq. (32) is due to Li^e' ) and dis¬ 
appears in the limit a —» oo. This implies a simple exponential 
decay of correlations as expected for a SR model. Conversely, 
(ii) the exponential decay is negligible when a < 1, as shown 
in the inset of Fig. 5(a), implying an essentially pure algebraic 
decay. 


D. Comparison with correlations of the LRI model 

In this section we compare the decay with distance of the 
correlation functions and ~L(i, j) of the LRK models 

with those of the correlations Cfj and of the LRI model, 
respectively, since they are related by a Jordan-Wigner trans¬ 
formation [see, e.g.. Sec. IIB2], The correlation functions 
Cfj and Cfj have been studied within the PM phase in 
Ref. [26] for a > 0.5. There, for a > 1 it was found that 
the long-distance behaviour is characterized by an algebraic 
decay with exponents y z — 2a and y x — a for the two 
correlations, respectively. 

For a > 1, our own calculations for CC and Cf are re¬ 
ported in Fig. 2(b) and Fig. 5(b), respectively. There, we 
compare the decay of correlations with that of the correspond¬ 
ing correlation functions in the LRK models, goiL j ) [see also 
Eq. (39)] and Z(), j), respectively, showing very good quan¬ 
titative agreement with the (semi-)analytic results. In agree¬ 
ment with Refs. [12, 53], an hybrid exponential and power- 
law behaviour is found for C“ [as well as for X(), j) and 
£ 20 , 7 )]. 

Conversely, when a < 1 only an algebraic tail is visible 
for the decay of Cf in Fig. 5 and C X C [as well as Z(z,))] 
in Fig. 2, since the initial exponential decay is too small 
to be observed, as expected from the discussion above. 
Moreover, no universal behaviour for the decay exponents is 
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identifiable for Cff and C zz . in this region. The exponents for 
C* X j and Si i. j) differ here in general, probably because the 
contribution of the string operators in //lri becomes more 
relevant. 

In the AM phase, instead, we find that the decay of C Z C 
[Fig. 6(a)] for a > 1 displays an algebraic tail with an expo¬ 
nent compatible with y z = a, mimicking Cff in the PM phase 
(a very precise estimate is forbidden by the decay oscillations 
between even and odd sites). These results are consistent with 
those of Ref. [12], obtained for the case a — 3. Notably the 
decay exponent is here always different from the value 2 a an¬ 
alytically computed for g 2 (i, j ) in Eq. (39). This discrepancy 
is again probably due to the role of the string operators, here 
quantitatively relevant even for a > 1. For a < 1 again no 
universal behaviour for y 7 is found. 


IV. EDGE MODES PROPERTIES 
A. Massless edge modes 





It is known that the SR Kitaev chain hosts modes localized 
at the edges [11] in the ordered phase for n/4 < 6 < 3n/4. 
These modes are fermionic and massless (in the limit L —> oo), 
thus they have Majorana nature [76] and are a consequence 
of the topological non-triviality of the ordered phase (see 
e.g. [77] and references therein). 

For the LRK models of Eq. (2), Majorana massless modes 
are found for a > 1 in the AMI region of the phase diagram 
of Fig. 1(b). Plotting the square of the wavefunction if/oj cor¬ 
responding to the zero edge gap Ao, defined in Sec. IIB 3, we 
find for e + 1 in Eq. (2) a hybrid exponential and algebraic 
decay with the distance from one edge of the chain. The ex¬ 
ponent of the algebraic decay of \if/oj\ 2 is found to be equal to 
2a. This behaviour is similar to that observed in Ref. [53] in 
the presence of LR pairing only. 

Interestingly, here the algebraic tail of the Majorana 


FIG. 7. (a) LRK model: spatial probability distribution for the Ma¬ 
jorana edge mode if/oj in the AMI phase for 6 - 0.357T, a = 3 and 
different e. When e = 1, the mode is exponentially localized at one 
edge, while a power-law tail gradually appears when e t 1. (b) LRI 
model: spatial probability distribution for the mode ipj 1 as function 
of the lattice site j in the AM phase for different a and 6 = 0.457T. 
The decay with j is numerically found to be exponential, (c) LRK 
model: Spatial probability distribution for the massive edge mode in 
the AM2 phase for 6 = O.ln for different a and e. The decay with 
distance is purely power law and the probability distribution is sym¬ 
metric with respect to half of the chain, (d) LRI model: Probability 
density pj = |w . \ + \wj \ of the two degenerate excited states 11 > 
and |2 - ) in the PM2 phase. When a > 1, Pj spreads into the bulk of 
the chain. When a < 1, Pj is exponentially localized near a single 
edge. In panel (b) and (d), we plot only the points compatible with 
the DRMG errors. 



FIG. 6. (a): Cff + , ;((|+R correlation for the LRI for 6 = 0.35tt (AM 
phase) and L = 200, computed from R 0 = 100. The data are fitted by 
the solid lines ~ l/R 7 -", y- reported in the plot, (b): C*L correlation 
at 6 = 0.207 n (PM phases) for a system of L - 30 sites. 


modes can be tuned to completely disappear by changing the 
parameter e that fixes the unbalance between the hopping and 
the pairing terms in the Hamiltonian of Eq. (2). Figure 7(a) 
shows an example for 6 = QAn and a - 3. When e + 1 the 
hybrid behaviour is fully visible and the state if/oj decays in 
the bulk with an algebraic tail. However, by approaching the 
value e = 1 the algebraic tail of (//(), decreases and eventually 
disappears. As a result, for e — 1 the wave function 
becomes exponentially localized at one edge. We find that 
this exponential localization for e — 1 is present also in the 
parameter region a < 1. Moreover, the edge gap Ao scales to 
zero exponentially with the system size L for all a as Fig. 8(c) 
shows. 

For the LR Ising Hamiltonian of Eq. (1), edge modes ap¬ 
pear in the AM phase for every a > 0. They have zero mass, 
since the edge gap £j 0 -) ~ Ljo-) [defined in Sec. IIB3] van¬ 
ishes. Examples of these edge modes for different a are given 
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A„ - LRK 




FIG. 8. (a) First 45 low energy states of the LRK models for L = 500 and 9 = O.ln as function of a. It is possible to distinguish the A 0 mode 
(solid points for three different e) separated from the A„>o bulk modes (light green triangles for e = 10). When a > 1 the vanishing of Ao 
for e t 1 signals the spontaneous breaking of the Z 2 symmetry of LRK and defines the AMI phase [see also the phase diagram in Fig. 1(b)]. 
When Ao > 0 the Z 2 symmetry is restored and the phase AM2 appears, (b) Main panel: A£(l) (empty circles) and AE(2) (solid triangles) 
for the LRI as functions of the inverse of the length of the chain L for 9 = 0.1597T and for three a reported in the plot. These two states are 
degenerate in the L —> oo limit, as finite-size scaling shows. Inset: A£(l) as function of a for 9 = 0.159tt. When a decreases, A£(l) increases, 
showing a change for the energy spectrum in the region PM2 of the phase diagram of Fig. 1(a). (c) Edge gap Ao for (squares) e = 0.5 and two 
different a as function of the system size L. When e = 1, the edge gap scales exponentially to zero for all a. For comparison, the case (crosses) 
e = 2 and a = 3 is also reported. There the edge mass displays an hybrid exponential and power-law decay. 


in Fig. 7(b), where we plot the square of the wavefunction ip ( P 
defined in Eq. (12). We find numerically that ip^ decays ex¬ 
ponentially with the distance from the edge of the chain for 
all values of a, as well as <p\ (not shown) from the opposite 
edge. 

B. Massive edge modes 

If we extend the analysis of the LR Kitaev Hamiltonian of 
Eq. (2) to different e and sufficiently small a, a totally new 
situation arises for the edge gap Ao and the edge modes: In 
the region denoted as AM2 in the phase diagram of Fig. 1(b), 
we find that Ay (which is zero for a > 1) becomes nonzero for 
a < 1 and e ± 1, also in the thermodynamical limit. 

This case is shown in Fig. 1(e), where we plot the edge gap 
Ao as function of a and 9 for e — 10. Between e = 1 and 
e - 10 we checked a continuous increase of the extension of 
the region AM2 and no transitions in between. Similarly in 
Fig. 8(a) we plot Ao together with several other single-particle 
energies as a function of a and for a fixed 9. These energies 
have been computed as described in Sec. IIB 3. In Fig. 8(a) 
the mass of the edge mode A 0 is easily recognizable for all e, 
since it is separated from all bulk modes by a finite gap. 

Consistently with the discussion above, for a > 1, Ao is 
zero as expected from the SR model, so that two degenerate 
ground states exist as the Z 2 symmetry of the model is spon¬ 
taneously broken. However, surprisingly, for a < 1 and e + 1, 
we find that Ao becomes finite and thus the ground state is 
unique. This indicates that the Z 2 symmetry of the model 
(which is broken for a > 1 in the AMI phase) is restored for 
cr < 1. As a consequence, the region AM2, where the Z 2 sym¬ 


metry is restored, must be separated from AMI by a quantum 
phase transition, even if no closure for the mass gap arises in 
the bulk. 

The wavefunction of the lowest massive Ao state is now 
given by the matrix element go; defined in Sec. IIB 3. By 
plotting the probability density |go;| 2 > we now obtain a 
localization on the edges that is symmetric with respect to 
the middle of the chain. This probability density decays 
algebraically when approaching half of the chain, as is clearly 
seen in Fig. 7(d). 

A similar wavefunction localization at the edges of the 
system is found also for the LRI model in the PM2 region 
in the phase diagram of Fig. 1(a) for a < 1, where the 
Z 2 symmetry is preserved. However, while for the LRK 
massive edge modes originate from Majorana edge modes 
present at a > 1, for the LRI the edge localization arises for 
excited states of the bulk spectrum in the region PM1. For 
a < 1, these states are degenerate, as shown in Fig. 8(b), and 
separated from the third excited state by a gap that is finite 
in the thermodynamic limit. Because of this degeneracy, 
we consider the probability density pj = \w . |" + |w. | , 
with vv. ’ wavefunctions of |(1,2)) defined in Sec. IIB3. A 
typical situation is depicted in Fig. 7(d), where we plot pj 
as function of the lattice site j for different values of a. For 
a > 1, pj is oscillating and delocalized in the bulk, while 
it is localized exponentially at the edges for a < 1 and is 
symmetric with respect to half of the chain. 

We leave as an open question whether the edge localization 
here signals the appearance of a new phase (and without mass- 
gap closure) with preserved Z 2 symmetry, similar to the LRK 
models above. 














13 


V. OBSERVABILITY IN CURRENT EXPERIMENTS 

Recent experiments with cold ions have made possible the 
realization of LR Ising-type Hamiltonians as Eq. (1) with 
0 < a < 3.5 [15, 16, 18, 20]. In these experiments, both static 
and dynamical spin-spin correlations, as well as the spectrum 
of quasi-particle excitations [25], can be measured with ex¬ 
treme precision, which in principle could allow for an analy¬ 
sis of some of the observables discussed above. For example, 
the observation of long-distance algebraic of correlations, as 
well as spectroscopic signatures of the formation of localized 
excited edge modes for a < 1 could allow for the precise de¬ 
termination of the properties of these LR models. 

One key aspect of experiments, however, is that experimen¬ 
tally attainable lengths for ion chains are currently limited to 
at most few tens of ions. It is thus natural to ask whether 
the characteristic long-distance decay of correlation functions 
described above can be observed in systems of such length. 
To explore this issue. Fig. 6 (right panel) shows the correla¬ 
tion C\ X R for Hamiltonian Eq. (1) in the PM1 and PM2 phases 
for a chain of L = 30 sites with open boundary conditions 
and for different a. For a > 2, the initial exponential decay 
dominates the correlations for R < 10, while a comparatively 
small algebraic tail is found for R > 10. For a < 1, however, 
the exponential part has essentially disappeared and the decay 
is purely algebraic at all distances, as expected from the dis¬ 
cussion of Secs. HID. This fundamental change of behaviour 
around a ^ 1 may be observable. We note, however, that the 
exponent y x of the algebraic decay is here different from that 
presented in Fig. 2(a), due to strong finite size effects in these 
systems. We find similar results for the correlation Cf R . 

On the other hand, the emergence of massive edge modes 
in the LRI chain could be a convenient diagnostic of the 
change of nature of the paramagnetic phase for a < 1. 


VI. SUMMARY AND OUTLOOK 

In this work we have analyzed the phase diagram of the 
long range anti-ferromagnetic Ising chain and of a class of 
fermionic Hamiltonian of the Kitaev type, with long-range 
pairing and hopping. We have clarified in what regions of 
the phase diagram violation of the area law occurs, and have 


provided numerical evidence and exact analytical results for 
the observed hybrid decay of correlation functions, which are 
found to decay exponentially at short range and algebraically 
at long range, for all a. We have further demonstrated the 
breaking of conformal symmetry along the critical lines in 
both models at low enough a. Most interestingly, for the 
fermionic models we have demonstrated for the first time 
that the topological edge modes can become massive for 
sufficiently small values of a < 1. This implies the existence 
of a transition to a novel phase without closure of the mass 
gap, to the known phase with massless Majorana modes for 
a > 1. We conjecture that the possibility of a phase transition 
with nonzero mass gap is due to the peculiar behaviour of 
LR correlations, showing power-law tails also when the gap 
does not vanishes. Similarly, we have found that excited 
bulk states in the paramagnetic phase of the Ising model can 
become localized at the edges of the chain for a < 1. 

This works may open several exciting research directions. 
The first question concerns the nature and topological prop¬ 
erties of the proposed new phase of the Kitaev model with 
a < 1, and of its localized edge modes. We conjecture that 
these massive edge modes are due to the hybridization of the 
Majorana modes at small a, due to the bulk overlap between 
their wave functions, whose decay is slower and slower for 
decreasing a. This aspect will be the subject of future studies. 
Another important open question is whether the appearance of 
massive edge modes may be connected also to the violation of 
the area-law for the entanglement entropy in these models. 

In general, these results represent counter-examples for the 
topological properties of existing topological models with 
long-range interactions, as recently analyzed in [50]. The 
question of whether a possible universal behaviour exists for 
topological models with long-range interactions is thus still 
wide open. 
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